%% Run this script to produce tables and figures of the EWM analysis
% This script does not involve the estimation of EWM results. 
% It only calls and summarizes previously estimated and saved point 
% estimates and confidence intervals to produce summary tables and
% treatment rule plots

clear; clc; 
cd('..\..')

savedir='.\intermediate_data\EWM_results';
savedir_tables='.\output\tables';
savedir_figures='.\output\figures';
savedir_results='.\output\matfiles';

addpath ('.\intermediate_data')
addpath('.\code\analyze')
addpath('.\intermediate_data\EWM_results')
addpath('.\output\figures')
%% Main specification
individual_ps=0; % =0 using treatment probability for each of the three waves
winter_prevuse=0; % =0 pre-treatment average calculated using specified number
% of months prior to treatment
pre_avg_n=12; % 12 month pre-treatment period
post_avg_n=12; % 12 month post-treatment period
wave=0; % =0 if pooled; otherwise =3, 6, or 7
bsimax=0; % number of bootstraps. =0 if no bootstrap

% Define tex output name
if (winter_prevuse)
    s_pre='pre_avg_winter_recent';
else
    switch pre_avg_n
        case 21
            s_pre='pre_avg_21';
        case 24
            s_pre='pre_avg_24';
        case 12
            s_pre='pre_avg_12';
    end
end

switch post_avg_n
    case 21
        s_post='post_avg_21';
    case 24
        s_post='post_avg_24';
    case 12
        s_post='post_avg_12';
    case 0
        s_post='all_post';
end
    
if individual_ps==1
    s_ps='ips';
else 
    s_ps='wps';
end

%fbase specifies output table prefix
fbase=sprintf('%s_%s',s_ps,s_pre,'_',s_post); 

if wave==0
    wave_name='pooled';
elseif wave==3
    wave_name='wave3';
elseif wave==6
    wave_name='wave6';
elseif wave==7
    wave_name='wave7';
end
%% Import sample 
sample_panel=readtable('elec_wave367_winterpre.csv');
% generate sample for analysis
[sample_cross] = generate_cross_sample(sample_panel,wave,winter_prevuse,pre_avg_n,individual_ps,post_avg_n);
clear sample_panel;

%% Produce PMC cost table 
cost = true; % switch: cost savings or energy conservation 
SMC=0; % use social marginal cost if SMC==1, use retail electricity price if SMC==0
tcost = 0.765; % cost of sending out letter for each household

ewm_pe_ci=cell(6,1);
covariates={'income' 'size' 'vintage' 'income' 'size' 'vintage'};
rules={'quadrant' 'quadrant' 'quadrant' 'cubic' ...
    'cubic' 'cubic'}; 
% need to manually specify based on bootstrap time stamp
s_now_batch_PMC={'20210306_142426' '20210306_152238' '20210306_171206' ...
    '20210312_204517' '20210314_083625' '20210315_090411'}; 
s_now_batch_delta_CI_PMC={'20210329_104706' '20210329_113026' '20210329_123937' ...
    '20210328_112141' '20210328_112933' '20210328_113650'};

[benchmark_1, ~, scaled_ATT, ~] = calculate_benchmark_savings(sample_cross,wave,tcost,SMC,bsimax);
benchmark_1_table = benchmark_1(1,:);
for i=1:6
    ewm_pe_ci{i}= generate_ewm_results_table(covariates{i},wave,winter_prevuse,SMC,fbase,scaled_ATT,s_now_batch_PMC{i},s_now_batch_delta_CI_PMC{i},rules{i},tcost);
end
savings_EWM = vertcat(ewm_pe_ci{:,:});
Rules= {'EWM-Quadrant';'';'EWM-Quadrant';'';...
    'EWM-Quadrant';'';'EWM-cubic';'';'EWM-cubic';'';'EWM-cubic';''};
savings_EWM= addvars(savings_EWM, Rules, 'Before', 'Covariates');
PMC_table=[benchmark_1_table; savings_EWM];

% export to latex
tablename=sprintf('PMC_wave%1.0f',wave);
filename=strcat(tablename,'_',fbase,'_fullCI','.tex');
save(fullfile(savedir_results,strcat(tablename,'.mat')),'PMC_table'); % save matfile 
table2latex(PMC_table, fullfile(savedir_tables,filename)) % export to latex


%% Produce kWh table 
cost = false; % switch: cost savings or energy conservation 
SMC=0; % use social marginal cost if SMC==1, use retail electricity price if SMC==0
s_cost = 'Kwh';
tcost = 0; % cost of sending out letter for each household

ewm_pe_ci=cell(6,1);
covariates={'income' 'size' 'vintage' 'income' 'size' 'vintage'};
rules={'quadrant' 'quadrant' 'quadrant' 'cubic' ...
    'cubic' 'cubic'}; 
% need to manually specify based on bootstrap time stamp
s_now_batch_kwh={'20210306_174428' '20210306_184312' '20210306_203213' ...
    '20210316_210520' '20210318_134122' '20220624_212327'}; 
s_now_batch_delta_CI_kwh={'20210330_044254' '20210330_052601' '20210330_063448' ...
    '20210328_115108' '20210328_115603' '20210328_120228'}; 

[benchmark_1, ~, scaled_ATT, ~] = calculate_benchmark_savings(sample_cross,wave,tcost,SMC,bsimax);
benchmark_1_table = benchmark_1(1,:);

for i=1:6
    ewm_pe_ci{i}= generate_ewm_results_table(covariates{i},wave,winter_prevuse,SMC,fbase,scaled_ATT,s_now_batch_kwh{i},s_now_batch_delta_CI_kwh{i},rules{i},tcost);
end
savings_EWM = vertcat(ewm_pe_ci{:,:});
Rules= {'EWM-Quadrant';'';'EWM-Quadrant';'';...
    'EWM-Quadrant';'';'EWM-cubic';'';'EWM-cubic';'';'EWM-cubic';''};
savings_EWM= addvars(savings_EWM, Rules, 'Before', 'Covariates');
kwh_table=[benchmark_1_table; savings_EWM];

% export to latex
tablename=sprintf('kwh_wave%1.0f',wave);
filename=strcat(tablename,'_',fbase,'_fullCI','.tex');
save(fullfile(savedir_results,strcat(tablename,'.mat')),'kwh_table'); 
table2latex(kwh_table, fullfile(savedir_tables,filename))


%% Produce SMC table 
cost = true; 
SMC=1; 
s_cost = 'SMC';
tcost = 0.765; 

ewm_pe_ci=cell(6,1);
covariates={'income' 'size' 'vintage' 'income' 'size' 'vintage'};
rules={'quadrant' 'quadrant' 'quadrant' 'cubic' ...
    'cubic' 'cubic'}; 
% need to manually specify based on bootstrap time stamp
s_now_batch_SMC={'20210306_210438' '20210306_220306' '20210306_235321' ...
    '20210321_083208' '20210323_082527' '20210324_113047'}; 
s_now_batch_delta_CI_SMC={'20210330_230517' '20210330_235355' '20210331_010925' ...
    '20210328_121412' '20210328_122040' '20210328_122709'}; 

[benchmark_1, benchmark_2, scaled_ATT, ~] = calculate_benchmark_savings(sample_cross,wave,tcost,SMC,bsimax);
benchmark_1_table = benchmark_1(1,:);

for i=1:6
    ewm_pe_ci{i}= generate_ewm_results_table(covariates{i},wave,winter_prevuse,SMC,fbase,scaled_ATT,s_now_batch_SMC{i},s_now_batch_delta_CI_SMC{i},rules{i},tcost);
end
savings_EWM = vertcat(ewm_pe_ci{:,:});
Rules= {'EWM-Quadrant';'';'EWM-Quadrant';'';...
    'EWM-Quadrant';'';'EWM-cubic';'';'EWM-cubic';'';'EWM-cubic';''};
savings_EWM= addvars(savings_EWM, Rules, 'Before', 'Covariates');
SMC_table=[benchmark_1_table; savings_EWM];

% export to latex
tablename=sprintf('smc_wave%1.0f',wave);
filename=strcat(tablename,'_',fbase,'_fullCI','.tex');
save(fullfile(savedir_results,strcat(tablename,'.mat')),'SMC_table'); 
table2latex(SMC_table, fullfile(savedir_tables,filename))

%% Plot EWM rules - Demographics 

% loop over three specifications: pmc, kwh, smc
cost={true false true}; 
SMC={0 0 1}; 
covariates = {'income','size','vintage'};
dimension="two-dim"; %two-dimensional rule or one-dimensional rule

for i=1:3
    if (cost{i})
        if (SMC{i})
            s_cost = 'SMC';
        else
            s_cost = 'PMC';
        end
    else
        s_cost='kwh';
    end

    if (cost{i})
        tcost = 0.765; 
    else
        tcost = 0;
    end
    
    for j=1:size(covariates,2) % loop over covariates
        covariate=covariates{j};
        % Plot treatment rules 
        [rules_sep_income, rules_comb_income, dens_map_income, heat_map_income] = summarize_savings_graph_covariate(covariate,wave,cost{i},winter_prevuse,SMC{i},dimension,fbase);
        % Save rule plot
        figname_rule=strcat('rule','_',fbase,'_',covariate,'_',s_cost,'_',wave_name,'.png');
        saveas(rules_comb_income,fullfile(savedir_figures,figname_rule))
        % Save heatmap
        figname_map=strcat('map','_',fbase,'_',covariate,'_',s_cost,'_',wave_name,'.png');
        saveas(dens_map_income,fullfile(savedir_figures,figname_map))
    end
end

%% Demographics rules for slides 
cost={true false true}; 
SMC={0 0 1}; 
covariates = {'income','size','vintage'};

for i=1:3
    if (cost{i})
        if (SMC{i})
            s_cost = 'SMC';
        else
            s_cost = 'PMC';
        end
    else
        s_cost='kwh';
    end

    if (cost{i})
        tcost = 0.765;
    else
        tcost = 0;
    end
    
    for j=1:size(covariates,2) % loop over covariates
        covariate=covariates{j};
        % generate rule plots without density 
        [rules_quad_income, rules_comb_income] = summarize_plots_for_slides(covariate,wave,cost{i},winter_prevuse,SMC{i},fbase);
        % quadrant rule only
        figname_rule=strcat('quad_rule','_',fbase,'_',covariate,'_',s_cost,'_',wave_name,'.png');
        saveas(rules_quad_income,fullfile(savedir_figures,figname_rule))
        % quadrant rule and cubic rule 
        figname_rule=strcat('rules','_',fbase,'_',covariate,'_',s_cost,'_',wave_name,'_','nodensity','.png');
        saveas(rules_comb_income,fullfile(savedir_figures,figname_rule))
    end
end

%% PMC, kWh, and SMC point estimates for ewm rules based on pre-treatment 
% consumption data only
cost={true false true}; 
SMC={0 0 1}; 
EWM_savings_baseline=cell(3,1);
EWM_univariate_savings=cell(3,1);
Scaled_ATT=cell(3,1);
scaled_ATT=cell(3,1);
ATE=cell(3,1);
 
 for i=1:3
      if (cost{i})
         if (SMC{i})
             s_cost = 'SMC';
         else
             s_cost = 'PMC';
         end
     else
         s_cost='kwh';
     end
 
     if (cost{i})
         tcost = 0.765; 
     else
         tcost = 0;
     end
     
     % Calculate benchmark savings in terms of ATE and ATT
     [Scaled_ATT{i},ATE{i}]=calculate_benchmark_savings(sample_cross,wave,cost{i},SMC{i},bsimax);
     scaled_ATT{i}=str2double(table2array(Scaled_ATT{i}(1,4)));
     
     covariates={'min' 'max' 'std' 'min' 'max' 'std'};
     rules={'quadrant' 'quadrant' 'quadrant' ...
    'cubic' 'cubic' 'cubic'}; 

    ewm_pe=cell(6,1);
    
    for j=1:6
    ewm_pe{j}= generate_ewm_results_table_pe_only(covariates{j},wave,winter_prevuse,SMC{i},fbase,scaled_ATT{i},rules{j},tcost);
    end
    
    [savings_EWM_one_dim_baseline] = generate_ewm_results_table_pe_only('income',wave,winter_prevuse,SMC{i},fbase,scaled_ATT{i},'onedim',tcost); 

    EWM_univariate_savings{i}=savings_EWM_one_dim_baseline;
    EWM_savings_baseline{i}=ewm_pe;
 end

 EWM_univariate_pe=vertcat(EWM_univariate_savings{:,:});
 pe_pmc = vertcat(EWM_savings_baseline{1}{:,:});
 pe_kwh = vertcat(EWM_savings_baseline{2}{:,:});
 pe_smc = vertcat(EWM_savings_baseline{3}{:,:});
 
 pmc_all=[Scaled_ATT{1}(1,2:5);EWM_univariate_pe(1,:);pe_pmc];
 kwh_all=[Scaled_ATT{2}(1,2:5);EWM_univariate_pe(2,:);pe_kwh];
 smc_all=[Scaled_ATT{3}(1,2:5);EWM_univariate_pe(3,:);pe_smc];
 
ewm_pe_all=nan(8,5);

ewm_pe_all=array2table(ewm_pe_all,'VariableNames',{'Rules','Covariates',...
    'kWh savings\ from\ EWM\ rules','PMC Savings\ from\ EWM\ rules','SMC savings\ from\ EWM\ rules'}); 

ewm_pe_all.Rules= {'Actual RCT';'EWM-Univariate';'EWM-Quadrant';'EWM-Quadrant';...
   'EWM-Quadrant';'EWM-cubic';'EWM-cubic';'EWM-cubic'};

ewm_pe_all.Covariates={'Scaled ATT';'mean';'min';'max';'std';'min';'max';'std'};
  
ewm_pe_all.(3)=kwh_all.(3);
ewm_pe_all.(4)=pmc_all.(3);
ewm_pe_all.(5)=smc_all.(3);

% EWM point estimates
tablename=sprintf('baseline_wave%1.0f',wave);
filename=strcat(tablename,'_',fbase,'_pe_all','.tex');
table2latex(ewm_pe_all, fullfile(savedir_tables,filename))

% Delta point estimates 
delta_ewm_rct_pe_all=nan(8,5);
 
delta_ewm_rct_pe_all=array2table(delta_ewm_rct_pe_all,'VariableNames',{'Rules','Covariates',...
    'PMC\ EWM\ v.\ RCT','kWh\ EWM\ v.\ RCT','SMC\ EWM\ v.\ RCT'}); 

delta_ewm_rct_pe_all.Rules= {'Actual RCT';'EWM-Univariate';'EWM-Quadrant';'EWM-Quadrant';...
   'EWM-Quadrant';'EWM-cubic';'EWM-cubic';'EWM-cubic'};

delta_ewm_rct_pe_all.Covariates={'Scaled ATT';'mean';'min';'max';'std';'min';'max';'std'};
  
delta_ewm_rct_pe_all.(3)=pmc_all.(4);
delta_ewm_rct_pe_all.(4)=kwh_all.(4);
delta_ewm_rct_pe_all.(5)=smc_all.(4);

delta_ewm_rct_pe_all=delta_ewm_rct_pe_all(2:8,:);

tablename=sprintf('delta_baseline_wave%1.0f',wave);
filename=strcat(tablename,'_',fbase,'_pe_all','.tex');
table2latex(delta_ewm_rct_pe_all, fullfile(savedir_tables,filename))


%% Plot EWM rules - Baseline consumption 
cost={true false true}; 
SMC={0 0 1}; 
covariates = {'min', 'max', 'std'};

for i=1:3
     if (cost{i})
        if (SMC{i})
            s_cost = 'SMC';
        else
            s_cost = 'PMC';
        end
    else
        s_cost='kwh';
    end

    if (cost{i})
        tcost = 0.765; 
    else
        tcost = 0;
    end

    % Plot treatment rules - one-dimensional 
    covariate="income";   
    dimension="one-dim"; 
    [rules_sep_one, rules_comb_one] = summarize_savings_graph_covariate(covariate,wave,cost{i},winter_prevuse,SMC{i},dimension,fbase);
    figname=strcat('rule','_',fbase,'_','one-dim','_',s_cost,'_',wave_name,'.png');
    saveas(rules_sep_one,fullfile(savedir_figures,figname))
    
    % Plot treatment rules - two-dimensional
    dimension="two-dim";
    for j = 1: size(covariates, 2)
        covariate = covariates{j};
        [rules_sep_min, rules_comb_min] = summarize_savings_graph_baseline(covariate,wave,cost{i},winter_prevuse,SMC{i},fbase);
        figname_rule=strcat('rule','_',fbase,'_',covariate,'_',s_cost,'_',wave_name,'.png');
        saveas(rules_comb_min,fullfile(savedir_figures,figname_rule))
    end
end

%% Summary table for delta ewm and rct (with CI) for baseline consumption specifications
cost={true false true}; 
SMC={0 0 1}; 
delta_univariate_table=cell(3,1);
delta_table=cell(3,1);
Scaled_ATT=cell(3,1);
scaled_ATT=cell(3,1);
ATE=cell(3,1);

s_now_baseline_all=cell(3,1);
s_now_batch_delta_CI_PMC={'20210329_130708' '20210329_191552' '20210330_013111' ...
    '20210328_113811' '20210328_113945' '20210328_114122'}; 
s_now_batch_delta_CI_KWH={'20210330_070235' '20210330_131325' '20210330_194325' ...
    '20210328_120358' '20210328_120527' '20210328_120709'}; 
s_now_batch_delta_CI_SMC={'20210331_014037' '20210331_075306' '20210331_141055' ...
    '20210328_122830' '20210328_122956' '20210328_123118'}; 

s_now_baseline_all{1,1}=s_now_batch_delta_CI_PMC;
s_now_baseline_all{2,1}=s_now_batch_delta_CI_KWH;
s_now_baseline_all{3,1}=s_now_batch_delta_CI_SMC;
s_now_univariate={'20210326_171103' '20210328_114428' '20210328_120816'};

 for i=1:3
      if (cost{i})
         if (SMC{i})
             s_cost = 'SMC';
         else
             s_cost = 'PMC';
         end
     else
         s_cost='kwh';
     end
 
     if (cost{i})
         tcost = 0.765; 
     else
         tcost = 0;
     end
     
     covariates={'min' 'max' 'std' 'min' 'max' 'std'};
     rules={'quadrant' 'quadrant' 'quadrant' ...
    'cubic' 'cubic' 'cubic'}; 

     % Calculate benchmark savings in terms of ATE and ATT
     [Scaled_ATT{i},ATE{i}]=calculate_benchmark_savings(sample_cross,wave,cost{i},SMC{i},bsimax);
     scaled_ATT{i}=str2double(table2array(Scaled_ATT{i}(1,4)));
     
    delta_pe_CI=cell(6,1);
    
    for j=1:6
    delta_pe_CI{j}= generate_delta_CI_table(covariates{j},wave,winter_prevuse,SMC{i},fbase,scaled_ATT{i},s_now_baseline_all{i,1}{1,j},rules{j},tcost);
    end
    
    [delta_univariate] = generate_delta_CI_table('income',wave,winter_prevuse,SMC{i},fbase,scaled_ATT{i},s_now_univariate{i},'onedim',tcost); 

    delta_univariate_table{i}=delta_univariate;
    delta_table{i}=delta_pe_CI;
 end


 delta_univariate_table=vertcat(delta_univariate_table{:,:});
 delta_table_pmc = vertcat(delta_table{1}{:,:});
 delta_table_kwh = vertcat(delta_table{2}{:,:});
 delta_table_smc = vertcat(delta_table{3}{:,:});
 
 delta_pmc_all=[delta_univariate_table(1:2,:);delta_table_pmc];
 delta_kwh_all=[delta_univariate_table(3:4,:);delta_table_kwh];
 delta_smc_all=[delta_univariate_table(5:6,:);delta_table_smc];
 
 Rules=["EWM-Univariate";"";"EWM-Quadrant";"";"EWM-Quadrant";"";...
   "EWM-Quadrant";"";"EWM-cubic";"";"EWM-cubic";"";"EWM-cubic";""];

delta_pmc_all=[Rules, delta_pmc_all];
delta_kwh_all=[Rules, delta_kwh_all];
delta_smc_all=[Rules, delta_smc_all];

tablename=sprintf('delta_baseline_wave%1.0f',wave);
filename=strcat(tablename,'_',fbase,'_pe_CI_pmc','.tex');
save(fullfile(savedir_results,strcat(tablename,'_PMC','.mat')),'delta_pmc_all'); 
table2latex(delta_pmc_all, fullfile(savedir_tables,filename))

tablename=sprintf('delta_baseline_wave%1.0f',wave);
filename=strcat(tablename,'_',fbase,'_pe_CI_kwh','.tex');
save(fullfile(savedir_results,strcat(tablename,'_kwh','.mat')),'delta_kwh_all'); 
table2latex(delta_kwh_all, fullfile(savedir_tables,filename))

tablename=sprintf('delta_baseline_wave%1.0f',wave);
filename=strcat(tablename,'_',fbase,'_pe_CI_smc','.tex');
save(fullfile(savedir_results,strcat(tablename,'_SMC','.mat')),'delta_smc_all'); 
table2latex(delta_smc_all, fullfile(savedir_tables,filename))

%% Summarize results on apply rules across waves (kWh, PMC, SMC) 
% only target wave 6 and wave 7
% need to run this block separately for coef=ewm and coef=delta
cost={true false true}; 
SMC={0 0 1}; 
savings_cubic_sample3=cell(3,1);
savings_cubic_sample6=cell(3,1);
savings_cubic_sample3_v2=cell(3,1);
savings_cubic_sample6_v2=cell(3,1);
scaled_att_sample_wave_weighted_cubic=cell(3,1);
scaled_att_target_wave_cubic=cell(3,1);
coef="delta"; % or coef="ewm"

for i=1:3
     if (cost{i})
        if (SMC{i})
            s_cost = 'SMC';
        else
            s_cost = 'PMC';
        end
    else
        s_cost='kwh';
    end

    if (cost{i})
        tcost = 0.765; 
    else
        tcost = 0;
    end
    
    [savings_cubic_sample3{i},savings_cubic_sample6{i},savings_cubic_sample3_v2{i},savings_cubic_sample6_v2{i},scaled_att_sample_wave_weighted_cubic{i},scaled_att_target_wave_cubic{i}] = ...
        generate_external_validity_table_cubic_all_specifications(sample_cross,winter_prevuse,s_cost,tcost,bsimax,SMC{i},savedir,fbase,coef);
end

% rules applied on original target wave
ev_cubic_table_target=[[[6;6;6],[3;3;3],round(savings_cubic_sample3{2,1},2),round(savings_cubic_sample3{1,1},2),round(savings_cubic_sample3{3,1},2)];...
         [[7;7;7],[6;6;6],round(savings_cubic_sample6{2,1},2),round(savings_cubic_sample6{1,1},2),round(savings_cubic_sample6{3,1},2)]];

ev_cubic_table_export_target=array2table(ev_cubic_table_target,'VariableNames',{'Target\ wave','Sample\ wave','Energy\ changes','Private\ cost\ changes',...
    'Social\ cost\ changes'});

if coef=="delta"
    filename_cubic=strcat('EV_cubic_table','_','delta','_',fbase,'_','target','.tex');
elseif coef=="ewm"
    filename_cubic=strcat('EV_cubic_table','_','ewm','_',fbase,'_','target','.tex');
end

table2latex(ev_cubic_table_export_target, fullfile(savedir_tables,filename_cubic))

% rules applied on sample wave * density ratio
ev_cubic_table_sample=[[[6;6;6],[3;3;3],round(savings_cubic_sample3_v2{2,1},2),round(savings_cubic_sample3_v2{1,1},2),round(savings_cubic_sample3_v2{3,1},2)];...
         [[7;7;7],[6;6;6],round(savings_cubic_sample6_v2{2,1},2),round(savings_cubic_sample6_v2{1,1},2),round(savings_cubic_sample6_v2{3,1},2)]];

ev_cubic_table_export_sample=array2table(ev_cubic_table_sample,'VariableNames',{'Target\ wave','Sample\ wave','Energy\ changes','Private\ cost\ changes',...
    'Social\ cost\ changes'});

if coef=="delta"
    filename_cubic=strcat('EV_cubic_table','_','delta','_',fbase,'_','sample','.tex');
elseif coef=="ewm"
    filename_cubic=strcat('EV_cubic_table','_','ewm','_',fbase,'_','sample','.tex');
end

table2latex(ev_cubic_table_export_sample, fullfile(savedir_tables,filename_cubic))

% scaled att sample wave weighted 
scaled_att_sample_wave_weighted_cubic_table=[[6;6;6],[3;3;3],cell2mat(scaled_att_sample_wave_weighted_cubic{2,1}(1,:))',cell2mat(scaled_att_sample_wave_weighted_cubic{1,1}(1,:))',cell2mat(scaled_att_sample_wave_weighted_cubic{3,1}(1,:))';...
    [7;7;7],[6;6;6],cell2mat(scaled_att_sample_wave_weighted_cubic{2,1}(2,:))',cell2mat(scaled_att_sample_wave_weighted_cubic{1,1}(2,:))',cell2mat(scaled_att_sample_wave_weighted_cubic{3,1}(2,:))'];
scaled_att_sample_wave_weighted_cubic_table=array2table(scaled_att_sample_wave_weighted_cubic_table,'VariableNames',{'Target\ wave','Sample\ wave','Energy\ changes','Private\ cost\ changes',...
    'Social\ cost\ changes'});
filename_att_cubic=strcat('EV_att_table','_','cubic','_','sample_wave','_','weighted','.tex');
table2latex(scaled_att_sample_wave_weighted_cubic_table, fullfile(savedir_tables,filename_att_cubic))

%% Summarize out-of-sample performance in pooled sample
pm=100;
reweight=0;
s_winter='';
sample_wave=0;
target_wave=0;
split_choices=[0.2,0.3,0.4,0.5];

cost={false true true}; 
SMC={0 0 1}; 
covariates={'income','size','vintage'}; 

% 50/50 split with 100 runs 
split=0.5;
% get summary tables
[perm_quadrant_summary_table, perm_cubic_summary_table] = generate_cross_validation_permutation_100_table(pm,split,reweight,covariates,cost,SMC,s_winter,target_wave,sample_wave,fbase,savedir_figures);

% combine summary tables
perm_summary_table = vertcat(perm_quadrant_summary_table, perm_cubic_summary_table);
rct_rows = find(string(perm_summary_table.Rule)=='RCT');
drop_row = rct_rows(2:end);
perm_summary_table(drop_row,:) = [];
s_split=sprintf('%d-%d',split*100,(1-split)*100);
tablename=sprintf('Average_ewm_cv_pooled_%ssplit_%dperms',s_split,pm);
filename=strcat(tablename,'.tex');
table2latex(perm_summary_table, fullfile(savedir_tables,filename))

%% Summarize EWM results with budget constraints 
constraint_type={'fixedshare';'capshare'}; 
fix_range=0;
covariates={'income', 'size', 'vintage'};

for i=1:size(constraint_type,1)
    constraint = constraint_type{i,1};
    constrained_table = summarize_budget_constraint_table(constraint,fix_range,covariates,cost,SMC,winter_prevuse,wave,fbase);

    switch constraint 
        case 'fixedshare'
            tablename=sprintf('Budget_constraint_%s_range_%0.2f',constraint,fix_range*100);
        case 'capshare'
            tablename=sprintf('Budget_constraint_%s',constraint);
    end        

    filename=strcat(tablename,'.tex');
    table2latex(constrained_table, fullfile(savedir_tables,filename))
end


%% Inequality analysis 
sample_panel=readtable('elec_wave367_winterpre.csv');
summarize_inequality_analysis(wave,sample_cross,sample_panel,savedir_figures,fbase)

